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1. Introduction 

The magneto-electric effect in insulating materials consists in the appearance of 
a finite, measurable electric bulk polarization P when a sample is subjected to 
an external magnetic field £?, and in the appearance of a finite, measurable bulk 
magnetization M when the sample is subjected to an external electric field E. The 
magneto-electric effect is quantified by the magneto-electric response tensor: 

° ij ~ dB, ~ dEi ' { ' 

where the derivatives are not necessarily taken at zero B or E fields. The effect 
has been observed in a variety of materials and its technological applications can be 
tremendous [21 [3]. 

Originally, the effect was sought in materials with broken time-reversal and 
inversion symmetries, but recently it was shown that certain classes of time-reversal 
or inversion invariant systems can display large magneto-electric responses [H [SJ [5] 
due to topological characteristics of their electronic structures. Such materials, 
which are now routinely synthesized and characterized in laboratories [3 [8], are 
called 3-dimensional topological insulators. It emerged in the recent years that 
the 3-dimensional topological insulators can be classified by their magneto-electric 
response. For example, for perfectly periodic, time-reversal invariant insulators, 
the magneto-electric response induces a Z 2 topological classification [4 and the 
topological part of the magneto-electric tensor, which can serve as the topological 
invariant for this classification, was shown [9] to equal the previously introduced 
[TUJ [TT], [T2] Z 2 invariants. This was an important development in the field of topological 
insulators because it provided a measurable bulk effect which sets apart the trivial and 
topological 3-dimensional time-reversal invariant insulators [13] . 

The magneto-electric response has ionic, spin and orbital contributions. Here 
we discuss only the orbital component, which, at its turn, can be decomposed into 
isotropic and traceless parts: 

1 3 

Olij = Sij - ^ a kk + &ij ■ (2) 
6 k=i 

Both components are experimentally measurable quantities but we will focus 
exclusively on the isotropic component: 

1 3 

a =~y2, ( xiii ( 3 ) 

6 3=1 

for it contains the topological response. For perfectly periodic, time-reversal or 
inversion invariant crystals at zero external fields, a was derived and explicitly 
computed in Refs. [H Q31 [5] for several lattice models of time-reversal or inversion 
symmetric insulators, showing that indeed a can take non-trivial quantized values. 
For perfectly periodic systems without time-reversal and inversion symmetries, 
Refs. [T3 [H] showed that a has additional contributions besides the topological 
component. 

Our goal here is to extend the above mentioned works to the case when disorder 
and magnetic fields are present, maintaining the assumption of a gap in the energy 
spectrum at the Fermi level. The main challenges arising in this more general setup 
are the breaking of the translational symmetry and the tedious linear response calculus 
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with respect to the magnetic field. Both challenges are overcome with the use of non- 
commutative geometric methods [T7] . The starting point for our work is Ref . [T] where 
non-commutative formulas for P and M have been reported. Given these results, we 
could have proceeded along two routs, one relying on the variation of P with B and 
the other one relying on the variation of M with E. We chose the first route because 
computing the variation of P with B fits better into the non-commutative geometry of 
the Brillouin torus |18| , which provides a set of explicit rules of calculus, including rules 
for how to compute the derivatives with respect to the magnetic field (see Sections [3] 
and [4|. These rules of calculus alone enabled us to complete the task at hand. 

For disordered time-reversal or inversion symmetric insulators (£? = 0), under the 
spectral gap assumption, we are able to formulate the following statement (similar to 
the existing one for perfectly periodic crystals [4j 5, 6 ): The total change Aapy] of a, 
induced by an adiabatic deformation along a path 7 in the model's parameter space, is 
equal to half of the second Chern number over the manifold (7 — £7) x Brillouin torus 
(when the units are chosen such that e = h = 1). Here, £ is either the time-reversal 
operation or the inversion operation and the path 7 starts and ends with a time- 
reversal or inversion symmetric system. Since Act [7] can change only by an integer 
number if the path 7 is re-routed, it follows that the integer or half-integer character 
of Aa[7] is independent of the deformation path and this extends the well-known Z 2 
topological classification to the disordered case. We would like to mention that the 
methods of non-commutative geometry enabled us to derive this result directly in the 
thermodynamic limit and to avoid the use of arguments based on twisted boundary 
conditions on large but finite volumes. The latter arguments are hard to control in the 
thermodynamic limit and they are problematic when strong disorder is considered. 

The explicit formula for a (Eq. 47) and its topological properties discussed in 
the present work is a small step in a program that aims at elucidating what happens 
when strong disorder is present. By strong disorder we mean situations when the 
Fermi level is embedded in a dense localized spectrum rather than being positioned 
in a spectral gap. The crucial question to answer is if the isotropic magneto-electric 
response remains quantized and invariant under the assumption of a mobility gap (as 
opposed to a spectral gap) at the Fermi level. A rigorous answer to this question 
will finally resolve the problem of existence (or lack) of robust extended states in 
3-dimensional topological insulators. 



2. Disordered lattice-models in the presence of magnetic fields 

The analysis reported in this work applies to generic disordered lattice-models in 
3-dimensions. The Hilbert space for such a model is spanned by vectors of the 
form |n, a), where n £ Z 3 is a node of the lattice and a — 1, ...,D labels the 
orbitals associated with that node. The scalar product (•, •) is by definition such that 
(n, a\m, (5) = S nrn S a /3. The generic lattice-Hamiltonians, in the absence of a magnetic 
field and disorder, take the form: 

^EE^-j^x™'/ 3 !' (4) 

n,a m,/3 

Accurate, material-specific lattice-models can be generated from ab-initio computa- 
tions using Wannier function representations or simply by fitting the measured band 
structure 19J . In almost all cases, the lattice-Hamiltonians involve only a finite num- 
ber of hoppings between one unit cell and its neighbors, which is something we will 
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assume in the followings. 

In the presence of disorder and a uniform magnetic field B, which will be 
parametrized as a 3 x 3 antisymmetric tensor B: 






B 3 


-B 2 


-B 3 







B 2 


-B x 






B = —B 3 B x | , (5) 
V B 2 - 
the lattice-Hamiltonian becomes: 

where the phase factor e 17r (" I - B "') encodes the effect of the magnetic field through 
the Peierls substitution (2U]. Above, (, ) denotes the 3-dimensional Euclidean scalar 
product and i = \J — f . The hopping amplitudes now have a random component, with 
a typical form: 

C-m( w ) = C-m + ^m.n, (7) 

where uJm,n are independent random variables uniformly distributed in the interval 
5]. The collection of all random variables lu = {uJm,n}m,n can be viewed as a 
point in an infinite dimensional configuration space S7, which can be equipped with a 
probability measure: 

dP(u) = Y[ du m , n . (8) 

in 11 

Furthermore, the discrete Z 3 additive group has a natural action on O given by: 

[Sa^)m,n ^m+a.ra+a? (9) 

which acts ergodically on SI and leaves dP(oj) invariant. If U a denotes the magnetic 
translation by a: 

U a \n,a) =e i7r{a ^ n) \n-a,a) 1 (10) 
then {-ff^l^gn forms a covariant family of operators in the sense that: 

UaH^Ua 1 = H iaL0 . (If) 

In the above settings, the quadruple (O, dP(cj), {t a } aS z 3 7 {H^Men) determines 
a homogenous system |18j . for which the intensive thermodynamic variables are, with 
probability one, independent of the disorder configurations. The results reported in 
this paper can be extended to more general models of homogeneous systems, but one 
should keep in mind the model of Eq. [6] as a prototypical example. 

We end this section with a comment about discrete versus continuum models. 
Definitely, working with lattice models reduces tremendously the technical complexity 
of the problem. One crucial difference is that the algebra of observables does have a 
unit in the discrete case, and it doesn't in the continuum case. Also, when working 
with discrete models, we are dealing only with bounded observable but this will not be 
the case for continuum models. Hence, extending the theory to the continuum models 
will require more work at the technical level (note also that the results of Ref. [T], 
which is starting point for our work, were established for discrete models). However, 
even if we do extend the theory to continuum models, an explicit computer calculation 
for a real material will not be possible, simply because of the sheer complexity of the 
problem. This is the main reason the community working on topological materials 
focuses almost entirely on lattice models. We should mention that, at least for the 
periodic case, the Wannier function technique generates lattice models which exactly 
reproduce the bundle of the occupied states of the continuum models. 



A Non-Commutative Formula for the Isotropic Magneto- Electric Response 



5 



3. The non-commutative Brillouin torus (Bellissard [18j) 

3.1. The C* -algebra of the covariant observables 

If for perfectly periodic systems one is primarily concerned with the translationally 
invariant operators, for homogenous systems one is primarily concerned with the 
covariant families of operators. All relevant response functions and thermodynamic 
functions for a homogeneous model can be computed within the algebra generated by 
the covariant families of operators. This algebra can be given a C*-algebra structure 
in the following way. 

For {-FLij-^go, a covariant family of operators, all the information is encoded in 
the matrix elements: 

f ai0 (Lj,n) = (Q,a\F U) \n,p), (12) 

because all the other matrix elements can be obtained from f a ^{uj,n) using 
translations. Hence, instead of working with covariant families of operators, one can 
work with functions: 

f :QxZ s ^Mdxd, (13) 

where A4 dxd is the space of D x D complex matrices, and with the algebraic 
operations that follows from the operator algebra when written for the matrix elements 
of Eq. [12J These operations are: 

(/ + 9)(u, n) = f(uj, n) + g(u), n), 

(14) 

(/ * g)(u>, n) = J2 mezd f(u>, m) 5 (t„>, n - m )e^ Bm ) . 

For the beginning, one should consider only continuous functions / with compact 
support. The algebra generated by such elements will be called Ao- One can easily 
go back from the algebra Aq to the algebra of bounded operators using the operator 
representation 7r w given by the formula: 

(W)l«)= E e^-^/fcVm-nJM, (15) 

where \n) and \m) should be viewed here as column matrices with entries |m, a) and 
|m,/3). The result of a representation {it u f}u)en is automatically a covariant family 
of operators. For example, the Hamiltonian of Eq. [6] is generated by the element: 

h a Au,n)=tf(w). (16) 

If one introduces the norm on Ao- 

||/|| = sup IK/H, (17) 
wen 

and the ^-operation: 

/*(o; ) n)=/(t-V-n)t ) (18) 



then the completion of _4o under the norm of Eq. 17 becomes a C*-algebra, which will 
be denoted by A. The C*-algebra structure allows one to develop the spectral problem 
and the functional calculus in a purely algebraic framework, without any reference to 
wavefunctions. The algebra A is the smallest algebra that one can consider and still 
be able to compute the response functions of the homogeneous systems. 
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3.2. The non-commutative differential calculus 

The classical differential manifolds arc completely determined by the algebra of 
smooth functions defined over the manifold and by the differential calculus with these 
functions. In the non-commutative setting, one replaces the commutative algebra of 
smooth functions over the classic Brillouin torus with the non-commutative algebra 
A, and then builds a non-commutative differential calculus over A. For the latter, one 
replaces fc-integration over the classical Brillouin torus by a trace over A: 

T(/)= f dP(e>) tr{f(w,0)}, (19) 

JS! 

where tr is the trace over Mdxd space, and the fc-derivations by linear automorphisms 
of algebra A: 

(djf)(u,n) =injf(u,n), j = l,...,d. (20) 

Then the triplet (A, T, d) determines a non-commutative manifold called the non- 
commutative Brillouin torus. 

The trace T has the standard properties: 

(i) Positivity: T(f * /*) > 0, 

(ii) Cyclic-city: T(f * g) = T{g * /). 

For the translational invariant case, T reduces to the ordinary fc-intcgration, and for 
the general case with disorder and a magnetic field, the physical meaning of T is that 
of the trace per unit volume: 

lim ^TtaIK/Xtt^) . . .} - T(f * g • • •), (21) 

A->oo |A| 

where Tta means the trace over the quantum states inside the box A. 

The derivations are not defined over the entire algebra A but only over the sub- 
algebra of "differentiable" functions: 

C\A) = {f G A, \\djf\\ < co, j = 1, 2, 3}. (22) 

More general, the sub-algebra of iV-th differentiable functions is defined as: 

C N (A) = {f € A, W^dZ'dg'fW < co, ai + a 2 + a 3 = N}. (23) 

We will only work with C°°(A) elements, so from now all elements will be assumed to 
be part of this sub-algebra. 

For convenience, we collect below a list of well-known rules of calculus (j, k = 
1,2,3): 

(i) The derivations commute: 

djd k f = dkdjf. (24) 

(ii) The derivations are ^-derivations: 

dj(f*) = (djf)*. (25) 

(iii) The derivations satisfy the Leibniz rule: 

9j{f*9) = {d j f)*g + f*{d J g). (26) 

(iv) The operator representations of the derivations are: 

n u} (d j f) = -i[x j ,7r u f], (27) 
where x — (x\, . . . , x s ) is the position operator. 
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(v) For / invertible in A, then: 

djf- 1 = -r 1 * (djf) * f~\ (28) 

(vi) If $ and <£' are any differentiable functions on the complex plane, then [21 : 

7W) * d^'(f)) = 0. (29) 

In particular: 

T(djf) = 0. (30) 
If combined with the Leibnitz rule, one gets the partial integration rule: 

T(d j f*g) = -T(f*d j g). (31) 

(vii) If $ is an analytic function in a neighborhood of the unit circle and u is any 
unitary element (u* = u^ 1 ), then [2"!2] : 

T{^{u)dju) = a- 1 T{u- 1 d j u), (32) 

where a_i is the coefficient of z^ 1 in the Laurent expansion = X)^L-oo a n,z n ■ 

4. The Ito derivative with respect to the magnetic field 

As we already mentioned, the traditional linear response calculus with respect to the 
magnetic field can be very difficult. In the context of the magneto-electric response, 
this was already highlighted in Ref. [23]. In the non-commutative setting, however, 
this task is quite straightforward, in great part due to the pioneering work by Rammal 
and Bcllissard [H]. According to this work, the derivative with respect to B of the 
expected values can be computed as: 

d B .T(f)=T(6 j f), j — 1, 2, 3, (33) 

where <5j's are the Ito derivatives, which can be computed using the following rules of 
calculus (with j taken mod 3): 

<V4/ ik*J. (34) 
djf* = (SjfT, (35) 

S j (f*g) = (S j f)*g + f*(S j g) (36) 
+ ^( d j+if * d j+2g - d ]+2 .f * d j+ ig), 

S 3 (r 1 ) = -.r 1 ^6 J f^f- 1 (37) 

These rules of calculus are enough to compute any response function involving the 
magnetic field. For example PQ: 

Proposition 1. If p is an idempotent element (p*p = p), such as the projector onto 
the occupied states for our model, then: 



and 



p * {5jp) *p = -^p* [dj+ip, d J+2 p] * p (38) 



(1 ~P) * (S 3 P) * (1 - j>) = ~(1 - p) * [d j+ ip, d j+2 p] * (1 - p). (39) 



A Non-Commutative Formula for the Isotropic Magneto- Electric Response 



Proof. According to the rule in Eq. [36] 

Sjp = Sj (p*p) = {5-jp) *p + p* (Sjp) + - [d j+ ip, d j+2 p] (40) 

Then Eq. 38 follows after one multiplies the above equation with p on both left and 
right sides and collects all three resulting p * Sjp * p terms. For the second identity, 
note that (1 — p) is also an idempotent, so we can replace p by (1 — p) in the first 
identity and use the fact that Sj(l — p) = —SjP- □ 

The off-diagonal components of Sjp do not have an expression that involves only 
p and its derivatives like in Eqs. [38] and [39] If p is the spectral projector of the 
Hamiltonian h onto the spectrum below the insulating gap: 

P=^<f(h-z)- 1 dz (41) 

where the contour encircles the occupied spectrum of h, then one will compute the 
off-diagonal part of 5jp using the rule in Eq. 37 on (h — z) . The result will depend on 
the particularities of h and will involve cross-gap matrix elements. In fact, in general, 
if we are dealing with a (analytic) function of the Hamiltonian, $>(h), such as the 
Fermi-Dirac function at finite temperatures, then we can use the representation: 

§(h) = ^-I$(z)(h-z)- 1 dz (42) 



and the rule in Eq. 37 on (h — z) , to compute 



Sj${h) = — jdz $(z) - z)- 1 * 5jh * (h - z)- 



+ 1 -(h-z)- 1 [d ]+1 h*(h-z)-\d 3+2 h*{h-z)- 1 }\. (43) 



5. The isotropic magneto-electric response 

The magneto-electric tensor is well defined only when referenced from a "standard" 
system. We will denote the reference system by ho and we will consider adiabatic 
time variations (i.e. infinitely slow) of the Hamiltonian h(t) which originate at this 
ho and terminate at some h. The main assumption is that the Fermi level is always 
in a spectral gap during the adiabatic variations. At the technical level, one needs 
the assumption that d t h(t) — at the beginning and at the end of the adiabatic 
deformation. Now, let p(t) be the projector onto the states below this spectral gap: 

p(t)= (l- X (h(t)-e F ))/2, (44) 

where x( x ) — sl S n ( x )- Under these conditions, Ref. pQ established the following result. 

Proposition 2. The change of the electric polarization during an adiabatic variation 
can be computed via the formula: 

AP j =ifdtr(p*[d t p,d j p]), (45) 

where the i-dependence of p is implied. 

It is worthwhile to add a few more lines about the work of Ref. [T] . For the context 
of disordered systems under magnetic fields, these results are important because 
they lift the (super-) adiabatic theorem [2S], including the full non-adiabatic error 
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estimation, to the non-commutative Brillouin torus. As such, one can foresee many 
additional applications which previously have been available only for the periodic case. 
Another notable point is that the formula for orbital magnetization derived in Ref. pQ 
does not require a spectral gap but only a mobility gap at the Fermi level. As such, 
there is hope that the theory of the magneto-electric response can be pushed into the 
regime of strong disorder. This is extremely important for the reasons to be discussed 
in the next sections. 

Starting from the result stated in Proposition 1, we shall prove the following 
statement. 

Proposition 3. In the above conditions, the variation in the isotropic magneto- 
electric response: 

= = dAP j/ dB i ( 46 ) 
3=1 

can be computed via : 

Aa = i J dt e a/3l sT(p * d a p * d p p * d 7 p * d s p) 



+ \r( X (h-e F )*d j p*5 j p) / (47) 



final 



initial 



3 

where summation over the repeating indices is assumed. The indices a, 5 run 
through t and j — 1, 2, 3, and e is the signature of the permutation. 



Remark. The first (second) term on the right hand side of Eq. 47 is called topological 
(non-topological) for reasons to be discussed shortly. When the magnetic field and 
disorder are turned off, the topological part of Aa reduces to the formula derived 
in Refs. [U [T51 [TB] in the context of perfectly periodic insulators. Showing a similar 
correspondence for the non-topological term is a more involved task and is deferred 
to future investigations. 

Proof. Since: 

3 

Aa = ijdt^J2 T & (P * d ?PV) > ( 48 ) 

3=1 

our main task is to compute: 

3 

£T(*i(P (49) 

3=1 

Using the rule of calculus from Eq. [36] we have: 

5j (p * [d t p, djp] ) = 8 j p* [dp, djp] +p*8 j [d t p, djp] 

+ 2dj+iP*dj+2[dtP,djp] - -d j+ 2p*dj + i[d t p,djp]. (50) 
The contributions to Aa from the last two terms above identically cancel out, as it 



can be easily seen by using an integration by parts (see Eq. 35 ). The Ito derivative of 
the commutator is: 

6j[d t p,djp] = [Sjdtp^djp] + [d t p,Sjdjp] 

+ \ {[dj+idtP, d j+2 djp] - [d j+2 d t p, d j+ idjp\). (51) 
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The contribution to the integrant in Eq. [48] from the first two terms of Eq. [51] is: 

3 
3=1 

3 3 

= 2 T([dtP, d jP ] *5 Pj ) + j Yl T( X (h - e F ) * d jP * S jP ) . (52) 
i=i 3=1 

Indeed: 

T(p * [Sj d t p, d jP ] +p* [d t p, Sj djp] ) (53) 
= T(p * Sjdtp * dj P — p * djp * Sjdtp + p * dtp * Sjdjp — p * Sjdj P * dtp) 
— T(djp * p * Sjdtp — p * djp * Sjdtp + p * dtp * Sjdjp — dtp * p * Sjdjp) 
= T(dt(dj P * p * djp) — dtdjp *p * djp — djp * dtp * Sjp 

- dt(p * djp * Sjp) + dtp * djp * 8j P + p * dtdjp * Sjp 

- dj (p * d t p) * S jP + dj (d tP *p)* Sjp) 

= d t T((dj P * p - p * djp) * 5jp) + 2T((d t p * d 3 p - d 3 p * d t p) * Sj P ) , 

and the affirmation follows because djp * p = (1 — p) * djp. 

The contribution to the integrant in Eq. [48] from the last two terms of Eq. [5l] is 

zero: 

3 

X] T (p * ([dj+idtp, d j+2 djp\ - [d j+2 d t p, dj+idjp})) = 0. (54) 

3=1 

Indeed, after an integration by parts, this contribution can be written as: 

3 

X T( - d-j+ip * [d t p, dj +2 djp\ -p* [d t , dj+^+^p] 
+ d j+2P * [d t p, d j+ idjp] +p*[d t , d j+2 d J+ idjp}) 

3 

= T( - d j+ ip * [d tP , d j+2 djp] + d j+2 p * [d t p, d ]+1 d jP \) 

3=1 

= T( - d 2 p * [d t p, d 3 d lP ] + d 3 p * [d t p, d^p] 

- d 3 p * [d t p, did 2 p] + dxp * [d tP , d 3 d 2 p] 

- dip * [d t ,d 2 d 3 p] + d 2 p * [d t ,did 3 p]), (55) 

and in the final lines one can easily identify pairs of terms that identically cancel out. 
So far the calculation stands as: 

3 3 

Y T ( 5 j(P * fttP, djp})) = 3 Y T([d t p, d jP ] * Sjp) 

3=1 3=1 
d 3 

3=1 



Now: 



T{[d t p,djp}*Sj P ) (57) 
T([d t p, djp] * Sjp * p) + T([d t p, d 3 p] * S d p * (1 - p)) 
T(p * [d t p, djp] * S^ *p) + T((l - p) * [d t p, djp] * Sjp * (1 - p)) , 
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where in the last line we used the cyclic-city of the trace and the idempotency of p 
and I— p. Given the fact that p and (1 —p) commute with [d t p, djp], we can continue 



as: 



. . . = T([%>, djp] *(p* Sjp *p+(l-p)* 5jp * (1 - p))) (58) 
at which point we can use the rules in Eqs. [38] and [39] to continue as: 

. . . = T([d t p, djp] *(--p* [d j+ ip, d J+2 p\ * p 

+ ^C 1 -P) * [ d 3 + lPi d 3+2P\ * (1 -P))) 

= ^T((2P - 1) * [dtp, d jP ] * [d j+1 p, d j+2 p}). (59) 
It remains to demonstrate that: 

3 

J2 T([d tP , d jP ] * [8 j+1 p, d J+2P \) = (60) 
which follows after an explicit display of the terms involved in the sum: 

3 

T{[d t p, djp] * [d j+ ip, dj+2p\) 

3=1 

= T[ptP * d\p * dip * dzp — dtp * dip * dsp * d 2 p 

— dip * dtp * d 2 p * dsp + dip * dtp * d^p * d 2 p 
+ d t p * d 2 p * d^p * d\p — d t p * d 2 p * d\p * d$p 

— d 2 p * dtp * d^p * dip + d 2 p * dtp * d\p * d%p 
+ dtp * d$p * dip * d 2 p — dtp * d$p * d 2 p * dip 

— dzp * d t p * dip * d 2 p + d^p * d t p * d 2 p * d\p) (61) 

and using the cyclic-city of the trace one can easily identify the pairs of terms that 
identically cancel out. At this point we arrived at: 

3 g 3 

X! T(Sj(p * [dtp, djp])) = - ^ T{p * [dtp, djp] * [d J+1 p, d J+2 p}) 

3 = 1 1 3=1 

+ ^^T( X {h-e P )*d jP *6 jP ), (62) 

3=1 

and the first term on the left side can be anti-symmetrized and brought to the desired 
form. This concludes the proof. □ 



6. Quantization and Invariance 

For the purpose of this discussion, it is more convenient to view the adiabatic time 
evolution from an initial to a final configuration as an adiabatic interpolation between 
two points in the parameter space of the model. The adiabatic interpolations will 
be represented as curves in this parameter space and the time-integral can be viewed 
as an integral along such curve. We automatically assume that the spectral gap at 
the Fermi level remains open during the adiabatic interpolations. We will use the 
notations 7 + 7' to denote the interpolation path obtained by joining the interpolation 
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Figure 1. (a) Various adiabatic interpolations between the initial (ho) and final 
(h) systems, (b) Here ho and h are time-reversal symmetric (TRS), but they 
are adiabatically interpolated via paths that are not necessarily time-reversal 
symmetric (©7 7^ 7). The shaded regions represents regions in the model's 
parameter space where the spectral gap closes. The TRS line is invariant to 
the time-reversal operation. 



paths 7 and 7', —7 to denote the interpolation path 7 when taken in reverse, and 
Aa[7] to indicate that the variation in the magneto-electric response depends on the 
interpolation path 7. 

The first question we want to address is what happens with the variation of the 
isotropic magneto-electric response when one considers different interpolation paths, 
7(i) and "f'(t), between the reference and the final systems. Since: 

Aa[7] - AafY] = Aa[7-7'], (63) 

we can write: 

Aa[-f] - Aa[y'] = ~<(,dt e aM ff{p * d a p * d p p * <9 7 p * d 5 p) , (64) 
7-7' 

and the important observation here is that the t-integral is now taken over a closed 
loop. One can then recognize on the right-hand side the expression of the second Chern 
number C2 over the manifold (7 — 7') x Brillouin torus, written in the real-space rather 
than in the fc-space. Because of the finite-gap assumption, the second Chern number is 
quantized, a fact that follows from a general argument by Connes [26] , which says that 
the pairing between the projectors and cyclic cocycles is a homotopy invariant. Hence 
the quantization occurring in the clean limit persists when the disorder is turned on, 
at least as long as the spectral gap remains open. 

Therefore, the answer to our first question is that the difference in Aa[7], when 
computed along different adiabatic evolutions, is always an integer number (we recall 
that our units were taken such that e = h = 1). We illustrate in Fig. 1(a) various 
possible situations. For example, if the interpolating paths 7(f) and Y(i) can be 
deformed into each other without closing the spectral gap, then Ci is necessarily 
zero and Aa[7] is identical for the two interpolation paths. But, if the interpolating 
paths cannot be deformed into each other without closing the spectral gap, such as it 
happens for j(t) and j"(t) in Fig. 1(a), then C2 may take non-trivial integer values so 
Aa[7] and Aa[7"] cannot be expected to be identical, but we can say with absolute 
confidence that the difference between the two will be an integer number. 

We now turn our attention to adiabatic interpolations between time-reversal 
invariant systems, in which case we need to turn the magnetic field off. Symmetry 
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considerations can be used to show that only the topological part of Aa survives for 
such such systems [151 116) . The time-reversal operation is given by the map: 



6 : .4| B=0 ^ -4| B=0 , (ef)(GJ,n)=e ms *f(u>,n)e-™\ (65) 

where s y is a matrix from Mdxd representing the y-component of the spin operator, 
and the line at the top represents complex conjugation. Some immediate properties 
of time-reversal map 8 are: 

e(/*<7) = (e/)*(e<?), (66) 

and 

djQf = -®(d 3 f). (67) 

Now, let Iiq and h be two time-reversal invariant Hamiltonians, Qho = ho and 
Qh = h. The key observation is that, if 7 (i) is an adiabatic interpolation between /i 
and h, then 07 is also an adiabatic interpolation between h$ and h (with a strictly 
positive spectral gap), and furthermore: 

Aa[ 7 ] = Aa[-9 7 ]. (68) 

This can be seen from the explicit formula of Eq. [47] and the properties of the time- 
reversal map stated in Eqs. [66] and [67] We already established that: 

Aa[ 7 ]- Aa[6 7 ] = C 2 [j - Qj], (69) 

and this together with Eq. |68| leads to an important conclusion for this section: 

Aa[ 7 ] = ic 2 [ 7 -e 7 ]. (70) 

This is a generalization of the Qi-Hughes-Zhang argument [3] to the case of disordered 
insulators. 

In Fig. 1(b) we illustrate several possible situations. If 7 can be deformed to a 
time-reversal invariant path, 67 = 7, then C2 = and Aa[ 7 ] is necessary zero. But 
if ho and h cannot be connected with a time-reversal invariant path without closing 
the spectral gap at the Fermi level, then 7 must be taken as in Fig. 1(b) so that to 
avoid the closing of the spectral gap. In this case, C 2 [ 7 — 67] can take non-trivial 
integer values and consequently Aa[ 7 ] can take non-trivial values, which are either 
integer or half-integer numbers. If we take a different interpolation path, such as 7' in 
Fig. 1(b), then Aa[ 7 '] may change but only by an integer number, as it was already 
established above. As such, for all possible adiabatic interpolation paths, Aa[ 7 ] is 
either an integer number or a half-integer number. In other words, the integer or 
half-integer character of Aa[ 7 ] is topologically stable. 



7. Conclusions 



The results of the last section can be further developed into a Z 2 classification of 
the 3-dimensional disordered time-reversal invariant insulators. For this, one can 
take the reference systems to be the insulators with zero hopping amplitudes between 
different lattice sites. These systems can be called "trivial" since they do not display 
any magneto-electric response as there are no Peierls factors when a magnetic field 
is turned on. As such, we can postulate a — for the trivial systems. Then, the 
time-reversal invariant insulators fall into two classes (here we assume that all time- 
reversal symmetric insulators can be adiabatically connected to a representative from 
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the trivial class): one for which a = integer and one for which a = half- integer, and we 
can say with absolute confidence that an insulator from the first class and an insulator 
from the second class cannot be connected by a time-reversal invariant interpolation 
path without closing the insulating gap. An insulator from a given symmetry class is 
generally called topological [6 if it cannot be adiabatically connected with its atomic 
limit (obtained by turning off the hopping terms between distinct lattice sites) via 
an interpolation path that respects the symmetry of the class and without closing 
the insulating gap. If we follow this definition, then the class of 3-dimcnsional 
disordered time-reversal invariant insulators with a = half-integer contains only 
topological insulators. Of course, our analysis does not rule out additional, more 
refined topological sub-classifications of either a = half or a = half-integer classes. 

A Z2 classification of the disordered inversion symmetric insulators is also possible 
via similar arguments. This is the case because the inversion operation satisfies the 
properties written in Eqs. [66] and [67] which are at the core of the Z2 classification of 
time-reversal insulators. However, for inversion symmetric insulators, it is explicitly 
known [S] that the Z2 classification is not the end of the story as it can be furtherly 
refined. 

We would like to point out that the computation of Section [5] remains valid even 
in the regime of strong disorder where the insulating gap is closed and is replaced by 
a mobility gap (provided we accept the formula Eq. 45 for the electric polarization). 
If a non-commutative theory of the second Chern number paralleling that of the first 
Chern number [21] can be developed, then the magneto-electric response can be used 
to topologically classify the time-reversal and inversion symmetric insulators in the 
strong disorder regime, which is one of the long-standing open problems in the field. 

The non-commutative formula Eq. [47] of the magneto-electric response can be 
efficiently evaluated on a computer by using the methods developed in Ref. [27) . 
As such, we now have a practical tool to explore, at least numerically, the physics 
of strongly disordered 3-dimensional topological insulators. The previous numerical 
attempts [28 to compute the 3-dimensional Z 2 invariants based on twisted boundary 
methods were limited to small systems due to the extreme complexity of such methods, 
even though these calculations involved a gauge invariant, highly efficient formulation 
of the Z2 invariant 29J. From our experience with other non-commutative formulas 
(such as those for the first Chern number and for the Kubo formula [301 EH EI] ) , we are 
quite optimistic that the non-commutative formula for magneto-electric polarization 
will enabled computations for much larger volumes where finer and more accurate 
results can be obtained. 
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